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but for other controllable DG technologies like gas turbines, it is not an easy job to provide a PDF to 
describe the generation schedule. On the other hand, the values of loads cannot be always described 
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the effectiveness of the proposed tool, it is applied to a realistic distribution system and the results are 
analyzed and discussed. 
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List of symbols 


SP, apparent forecasted value of peak load in bus i 

Sy apparent fuzzy value of demand in bus i and 
demand level h 

PE active power of a DG unit in bus i in demand level h 

Pe iq active power imported from grid in demand level h 

so apparent power magnitude passing through sub- 
station in demand level h 

Yi admittance magnitude between buses i and j 

Oi admittance angle between buses i and j 

Pax capacity limit of existing feeder l 

c2 capacity of DG unit installed in bus i 

ST ax capacity limit of existing transformer 

I current magnitude passing through feeder l in 


demand level h 


Vet cut out speed of the wind turbine 

v, cut in speed of the wind turbine 

DLFp demand level factor in demand level h 

pwind generated power of wind turbine in bus i and 

demand level h 

Vinin minimum operating limit of voltage 

Vraz maximum operating limit of voltage 

Pat net active power injected to bus iin demand level h 
ne net reactive power injected to bus iin demand level 

h 

Np number of buses in the network 

Ne number of feeders in the network 

Nn number of demand levels 

OF reactive power of DG unit in bus iin demand level h 

Vrated rated speed of the wind turbine 

pwina rated power of wind turbine installed in bus i, 

k shape factor of weibul PDF of wind speed 

c scale factor of the weibul PDF of wind speed 

Vij voltage magnitude of bus j in demand level h 

hj voltage angle of bus j in demand level h 


1. Introduction 


Distributed generations (DGs) are defined as electric resources 
interconnected to the distribution networks [1]. The integration of 
DG units in distribution network was highly studied in the last 
decade. There are many technical [2], socio-economic [3] and 
environmental [4] reasons behind this trend. In this paper, the focus 
is on the effect of DG units on active loss of distribution networks. 
The effect of conventional DG technologies on active power losses is 
investigated in the literature. In [5], an efficient method is proposed 
for the determination of optimal sizes and location of DG units by 
using the Kalman filter algorithm. In [6], a goodness factor is 
proposed which implies the contributions of a DG unit to active and 
reactive power losses in the distribution system are modeled using. 
In [7], a MINLP-based optimization model is proposed to determine 
the appropriate location of DG units in a distribution network in 
order to decrease active losses and fuel costs. In [8], the optimal 
location and size of DG units are obtained using a highly efficient 
ordinal optimization method. In these methods, the considered DG 
technology is controllable and non-renewable. The technological 
development and the importance of using clean technologies have 
made the renewable energies more fascinating for distribution 
network operators (DNOs) specifically because they are inexhaust- 
ible and nonpolluting. These technologies include hydro, wind [9], 


solar [10], biomass [11] and tidal. Among these renewable energies, 
wind technology has evolved very rapidly over the past decade and 
reduction of capital costs, improvement of reliability, and efficiency 
have made the wind power able to compete with conventional 
power generation [12]. The renewable DG technologies like wind 
have special characteristics due to their main source of energy. 
Obviously, the primary energy source of a wind turbine is wind. The 
wind speed is not a constant quantity during the operation of wind 
turbine and is highly dependent on climate condition of the area 
which wind turbine is installed there. It is the reason that they 
exhibit uncertainty and variability in their output [13]. In the 
literature, some methods are proposed to model the impact of these 
uncertainties on distribution network performance. In[14], different 
scenarios are constructed based on the PDF of uncertain values and 
then a method is proposed to determine the optimal combination of 
different renewable technologies for minimizing active losses. In 
[15], a powerful tool was proposed based on Monte Carlo Simulation 
for modeling the uncertainties in the locations, exported energy and 
penetration level, the states (on/off) of the DG units. It is assumed 
that the mentioned uncertainties follow a probability distribution 
function and the mentioned PDF is available for DNO. In some cases, 
there is no PDF available for describing the behavior of an uncertain 
parameter. In these cases, the Fuzzy arithmetic [16] can be used. In 
[17,18], a fuzzy power flow model is proposed for modeling the 
uncertainty of loads ina network. They use the extension principle of 
fuzzy to find the membership function of output when the input 
values are fuzzy and their membership functions are available. The 
contribution of this paper is as follows: 


e An evaluation tool is proposed for DNOs which helps them to 
deal with the effect of stochastic (probabilistic) and fuzzy 
(possibilistic) uncertainties of renewable and conventional DG 
technologies on active energy losses, simultaneously. 


To do this, the introduced principle of fuzzy load flow in [18] is 
used to model the uncertainties of loads and generation of non- 
stochastic DG technologies. For modeling the stochastic behavior 
of wind turbines, the Monte Carlo Simulation is used. 

The paper is structured as follows: Section 2, describes the 
possibilistic method used for uncertainty modeling, Section 3, 
presents the concepts of Monte Carlo Simulation for probabilistic 
uncertainty modeling; the problem formulation is described in 
Section 4. The proposed solution algorithm is presented in Section 
5. In Section 6, the proposed model is applied on an actual 
distribution network and the simulation results are given and 
discussed. Finally, conclusions are drawn in Section 7. 


2. Possibilistic uncertainty modeling 


In possibilistic evaluation frameworks, for each uncertain value, 
i.e. Â, a membership function, i.e. u4(x), is defined which describes 
that how much each element, i.e. x, of universe of discourse, i.e. U, 
belongs to A. Different types of membership functions can be used 
for describing the uncertain values. Here, fuzzy trapezoidal 
numbers (FIN) with a notation A = (amin, 4L, Au, Gmax) are used 
as shown in Fig. 1. 


2.1. a-Cut method 


In engineering problems, the evaluation of a certain quantity is 
usually in the form of a multivariable function namely, y = f(x1, ..., 
Xn), if & are uncertain then y will be also uncertain, 
Y= f(X,...,Xn). The question is that, knowing the membership 
functions of uncertain input variables X%,, what would be the 
membership function of y. The a-cut method [16] answers this 
question in this way: for a given fuzzy set A defined on universe of 
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Fig. 1. Fuzzy trapezoidal number. 


discourse, i.e. U, the crisp set A” is defined as all elements of U 
which have membership degree to A, greater than or equal to a, as 
calculated in (1): 

AX = {xe Ulua (X) >a} a) 
A% = (A*,A*) 
The a-cut of each input variable, i.e. x¥, is calculated using (1), then 
the a-cut of y, i.e. y“ is calculated as follows: 


y® = min f(X%) (2) 


X“ € (X%_X*) (3) 


This means for each a-cut, two optimizations are done. One 
maximization for obtaining the upper bound of y“, i.e. y“, and one 
minimization for obtaining the lower bound of y“, i.e. y“. 


2.2. Defuzzification 


The defuzzification is a mathematical process for converting a 
fuzzy number into a crisp one [16]. In this paper, the centroid method 
[19] is used for defuzzification of fuzzy numbers. The defiuzzified 
value of a given fuzzy quantity, i.e. A, is calculated as follows: 


«fal xdx 
A= [ede (4) 


3. Probabilistic uncertainty modeling 


Some of the uncertain input parameters follow a of probability 
distribution function (PDF), such as the value of wind which follows 
a weibull PDF [20]. Monte Carlo Simulation (MCS) is a powerful tool 
for analyzing the uncertainties which follow any PDF. 


3.1. Monte Carlo Simulation 


The main concept of MCS method is described as follows: 
suppose a multi-variable function, namely y, y = Z), where Z = [Z,, 
. +) Zm], in which z, to Zm are random variables with their own PDF. 
The problem is, knowing the PDFs of all input variables, i.e. Z1 to Zm, 
what would be the PDF of y? The MCS acts as follows [21]: first of all, 
it will generate a value, i.e. z?, for each input variable z; using its own 
PDF and form the Z° = [z‘,... ,Z%,] and then calculates the value of y° 
using y°=f(Z°). This process will be repeated for a number of 
iterations. The trend of the output, i.e. y, will determine its PDF. 


4. Problem formulation 


The assumptions for modeling the two types of uncertainties, 
constraints and the objective functions are described as follows. 


4.1. Uncertainty modeling 


As already explained, the uncertain parameters are divided into 
two groups: possibilistic and probabilistic. In possibilistic uncer- 
tainty group the value load in each bus, DG generation which are 
not stochastic (controllable with decisions of their owners). The 
second group contains the stochastic generation of wind turbines 
which is probabilistically modeled. The description of the 
parameters of each group is as follows: 

Possibilistic parameters: 


e Load: It is assumed that the DNO can just describe it with a 
membership function: 


Shi = Spy x DLF, x (Dmin; DL, Du, Dmax) (5) 


where S?, is the apparent forecasted value of peak load in bus i and 
DLF, is the demand level factor in demand level h which takes 
values between 0 and 1. Finally, SP; is the fuzzy value of demand in 
bus i and demand level h. 


e DG generation pattern: The amount of energy which a controlla- 
ble DG unit injects into the network is uncertain and usually it 
depends on the decisions of DG owner so the DNO cannot have a 
PDF of it if there is not much historic data about it. The output 
power of a controllable DG unit is modeled using a membership 
function as follows: 


Pe h Caini E C0 Sia (6) 


where Cc is the capacity of DG unit installed in bus i and Be is the 
active power of a DG unit in bus i in demand level h. 
Probabilistic parameters: 


e Wind Turbine generation pattern: The generation schedule of a 
wind turbine highly depends on the wind speed in the site. The 
variation of wind speed, i.e. v, can be modeled using a weibull 
[22] PDF and its characteristic function which relates the wind 
speed and the output of a wind turbine [23]: 


wu- (O (0) m 


where k is the shape factor and c is the scale factor of the weibul 
PDF of wind speed in the zone under study. 

The generated power of the wind turbine is determined using 
its characteristics as follows: 


0 if v < v orv > ve 
i V—Vin pwind ; 
pe > rey as if Vratea < V < Vout (8) 
out in 
pind else 


where ping is the rated power of wind turbine installed in bus i, 
pwind is the generated power of wind turbine in bus i and demand 
level h, vs, is the cut out speed, v{, is the cut in speed and Vyareq is 


the rated speed of the wind turbine. 
4.2. Constraints 


The constraints that should be satisfied are described in this 
section, as follows. 
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4.2.1. Power flow constraints 
The power flow equations must be satisfied in each demand 
level dl and year t, is as follows: 


net D wind 
Phi -Pp + SEPE ` Phi 
dg=1 wind=1 


a = -aP + So a 
al 
Np 
Pry = Vaid YuVny cos(bni — ôn j — ij) 


o= Vh, DDAA sin(ôn i = ôn j — 8ij) 
j=1 


where Pris‘ and Q}$' are the net active and reactive power injected 


to the network in bus i. 


4.2.2. Operating limits of DG units 
The power factor of DG unit is kept constant in all demand 
levels as follows: 


Pe 


TECI 


cos o% = 


= const (10) 


4.2.3. Voltage profile 
The voltage of each bus in each demand level h should be kept 
within the safe operating limits: 


V min < Vii < V max (11) 


Vmin and Vmax are the minimum and maximum permissible limits 
of voltage, respectively. 


4.2.4. Thermal limits of feeders and substation 

To maintain the security of the feeders and substation, the flow 
of current passing through them should be kept below their 
thermal limit in each demand level. This constraint is described as 
follows: 


Th < lax (12) 


grid tr 
S S Simak 


where It is the current magnitude of feeder l in demand level h, If, 
is the thermal limit of feeder 4, ae | is the apparent power imported 
from main grid in demand level h and finally Sf., is the thermal 
capacity of substation. 


4.3. Active losses 


The total active loss of the network is equal to the sum of all 
active power injected to each bus, as follows: 


Nn / Np . 
Loss = S~ (Se + me) x Th 


h=1 \i=1 

Loss € (loSS in, loss,, lossy, loSSmax) 

a=0 

losSmin = Min Loss 

losSmax = Max Loss (13) 
a=1 

loss, = min Loss 

lossy = max Loss 

Subject to : 

(5) to (12) 


where Ty, is the duration of demand level h and N; is the number of 
demand levels. 


5. Proposed model for mixed possibilistic-probabilistic 
uncertainty modeling 


The DNO is always faced with situations where the nature of 
the uncertainties are neither pure probabilistic nor possibilistic. 
In such cases, a tool is needed to handle the uncertainties. In this 
paper an evaluation tool is proposed to deal with such situations 
as follows: the DNO has an index function which is multivariable, 
i.e. f(X, Z), where the possibilistic uncertain parameters are 
represented by vector X and probabilistic uncertain values are 
given by vector Z. To deal with such variables they are 
decomposed into two groups and are dealt with separately as 
the following steps: 


e Step 1: For each z; € Z, generate a value using its PDF, i.e. z? 
e Step 2: Calculate (y*)° and (y*)* as follows: 


(y* )° = min f(Z°,X%) 
(y*)° = max sie X”) (14) 


These steps are repeated to obtain the PDF of the parameters of the 
outputs memebership function. In fact, when both types of the 
uncertainties exist in the input variables then the parameters of 
the membership function are stochastic. 


e 
$ 
4 
Cady, 
r 
. F. 


Fig. 2. Geographical view of the actual distribution system under study. 
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Table 1 

Data used in the study. 
Parameter Unit Value 
c 8.78 
k 1.75 
Vmin Pu 0.95 
Vmax Pu 1.05 
Via m/s 3 
Vrated m/s 13 
Vou m/s 25 
Duin 0.850 
Dı 0.925 
Du 1.075 
Dias 1.150 
Smin 0 
Ss. 0.9 
Su 1 
Smax 1 
Nn 24 

Table 2 


DG capacities and locations. 


DG technology Bus (i) DG capacity (MVA) 
Gas turbine 15 0.5 

163 0.5 

283 1 

495 3.5 
Wind turbine 426 0.5 

344 0.5 


6. Simulation results 


The proposed methodology is applied to a 574-bus realistic 
distribution network which is shown in Fig. 2. This network 
consists of a 20kV substation with capacity limit, i.e. Stax = 
20 MVA and 573 feeders with 180 load points. 

The average load and power factor at each load point are 
55.5 kW and 0.9285, respectively. There are six DG units present in 
the network where four of them are dispatchable (by non-DNO 
entities) and two of them are wind turbines. The shape and scale 


F 
2 
o 


Demand Level Factors (DL 
o 
oa 


101112131415 161718192021 222324 
Hours 


= 
NL 
wL 
AL 
OL 
OL 
NE 
ool. 
oL 


Fig. 3. The variations of DLF, in each demand level. 


factors of the weibull PDF of wind speed and the other simulation 
parameters are given in Table 1. The DG capacities and their 
location in the network are given in Table 2. 

It is assumed that there are 24 demand levels in each year with 
equal duration of Tp = 365 h. The variations of demand level factors 
are depicted in Fig. 3. 

The proposed algorithm is applied to the introduced network 
and the total active losses are obtained. The total active losses of 
the network can be presented in two ways: the first method is 
representing it by a trapezoidal fuzzy number in which each 
parameter of the membership function is a probabilistic quantity 
and has a PDF or histogram. In Fig. 4, the probability distribution 
functions of four parameters are shown. These parameters, i.e. 
(l0SSmin, loss;, lossy, losSmax), describe the total loss in the 
network. 

The histogram of variation of each parameter is shown in Fig. 5, 
describes the distribution of samples in the total Monte Carlo 
experiments which has been 20,000 experiments in this study. 

The second method for representing the total active loss is 
calculating the crisp value of active loss using (4) in each Monte 


3400 3450 3500 3550 


PDF(x) 


lossy (MWh) 


4400 4500 4600 4700 
loss, ox (MWh) 


Fig. 4. Probability distribution functions of total loss membership function. 
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Fig. 5. Histogram of total loss membership function. 


Carlo experiment and then obtaining the distribution of this 
quantity as shown in Fig. 6. 

The cumulative distribution function of crisp value of total 
losses is also shown in Fig. 7. 

For example if DNO wants to know the probability of having 
more than a specific loss in their network. For example, that might 
be a question that what is the probability of having more than 
3100 MWh loss in the network? Referring to Fig. 7, the probability 
that total loss exceeds 3100 MWh, is equal to 0.4388. Further 
statistical information of total active losses, is given in Table 3, 
including the minimum, maximum, mean value and the standard 
deviation of them. 

The developed hybrid tool attempts to overcome limitations in 
evaluating the network losses when different sources of uncertainty 
exist. It is specifically applicable when renewable and conventional 
DG units are present in the network. This hybrid approach makes the 
DNO enable to evaluate active losses when there are stochastic 
generations and also controllable DG units in the network. The 


3080 3090 3100 3110 3120 3130 3140 3150 
Crisp value of Loss (MWh) 


Fig. 6. Histogram of crisp value of total losses. 


simulation results can help the DNO to have estimation about the 
amount of money he should pay for compensation of active losses. 
Although the proposed analysis is offline and the running time of the 
algorithm is not of major concern but the computation burden of the 
proposed algorithm can be reduced using load flow techniques 
developed for radial distribution networks and also the variance 
reduction methods proposed for reducing the number of necessary 
Monte Carlo experiments. This will highly increase the speed of the 
proposed algorithm. 


Crisp Values of Loss 


Comulative Probability Distribution 
Function of Loss 


0 f i i i i i 
2800 3000 3200 3400 3600 3800 4000 4200 
Crisp value of loss (MWh) 


Fig. 7. Cumulative distribution function of crisp value of total losses. 


Table 3 

Statistics of values of total active loss in MWh. 
Quantity Min Max Mean Std 
losSinin 301.23 3772.9 1046.2 68.748 
loss, 2540.5 3911.9 3481.2 50.169 
lossy 3626.3 4477 3884.2 20.209 
losSmax 3976.4 4677.7 4477.9 68.704 
Crisp loss 2837.1 4060.3 3101 34.91 
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7. Conclusion 


A method combining probabilistic and possibilistic uncertainty 
modeling is proposed for evaluation of active losses in the 
distribution network. In doing so, it uses the Fuzzy (a-cut) and 
Monte Carlo Simulation to model the uncertainties. Its use would 
be to enable DNOs to evaluate the effect of different DG 
technologies on the technical performance of the distribution 
network. The model considers probabilistic presentation of wind 
speed using a weibull PDF and possibilistic description of loads and 
generation pattern of dispatchable DG units for modeling the 
uncertainties. The proposed tool is implemented on a realistic 
distribution network to demonstrate its ability. 
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